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Abstract 

A high-performance Fortran code is developed to calculate the spin- and parity-dependent shell model nuclear level 
densities. The algorithm is based on the extension of methods of statistical spectroscopy and implies exact calculation 
of the first and second Hamiltonian moments for different configurations at fixed spin and parity. The proton-neutron 
formalism is used. We have applied the method for calculating the level densities for a set of nuclei in the sd-, 
pf-, and pf + gg/2 - model spaces. Examples of the calculations for 28 Si (in the sd-model space) and 64 Ge (in the 
Pf + g9/2-model space) are presented. To illustrate the power of the method we estimate the ground state energy of 
64 Ge in the larger model space pf + gg/2, which is not accessible to direct shell model diagonalization due to the 
prohibitively large dimension, by comparing with the nuclear level densities at low excitation energy calculated in the 
smaller model space p f. 
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1. Introduction 

The nuclear reaction theory requires exact knowledge 
of nuclear level densities. In the majority of cases, espe- 
cially those important for reactions at stellar conditions, 
experimental information on the excited nuclear states 
is not sufficient and the calculations are typically using 
the estimates based on the Hauser-Feshbach approach 
01, where the level density for specific quantum num- 
bers J" of nuclear spin and parity is a necessary ingredi- 
ent. The reaction rates can be very sensitive to the level 
density, especially in the Gamow window of excitation 
energies around the particle threshold ||2|[3]]. 

The theory of the nuclear level density has a long 
history starting with the combinatorial calculation by 
Bethe [4] that leads to the the back-shifted Fermi gas 
approximation [5, 6] that was improved over the years. 
The modern approaches going essentially in the same 
direction [7-9] use, instead of the Fermi-gas, an inde- 
pendent particle model in a mean field. It is understood 
that many-body correlations may significantly change 
the resulting level density, in particular moving down 
certain families of vibrational and rotational states. At 



the same time, pure collective models used on top of 
the single-particle excitations may suffer from double 
counting and have to be cut off in some way. 

The many-body correlations are fully accounted for 
by the exact diagonalization of the shell-model Hamil- 
tonian. Here the problems come from the prohibitively 
large dimensions that makes necessary to truncate the 
orbital space, a step not always well controlled. Alter- 
natively, one can use Monte-Carlo techniques [10 -15] , 
or other methods of statistical spectroscopy ifm 171. 
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Most of these methods [2,10-12,18,19] calculate the to- 
tal density of states and later use a spin-weight factor 
that includes an energy-dependent cut-off parameter to 
extract the level density for specific quantum numbers 
of spin, parity and isospin. Although there are recent ef- 
forts to improve the accuracy of such parameterizations 
II 1 911 - it was shown that the cut-off parameter has very 
large fluctuations at low excitation energy 112011 . The par- 
ity is usually taken as equally distributed, although there 
are attempts ifTil I2H1 to model the effect of the uneven 
parity-dependence of the level densities at excitation en- 
ergies of interest for nuclear astrophysics. 

Recently, we developed a consistent approach [20,22- 
25] to calculate the spin- and parity-dependent shell- 
model level density. The new effective interactions 
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for the appropriate model spaces are developed starting 
with the G-matrix J26ll and fixing the monopole terms 
or/and linear combinations of two-body matrix elements 
to experimental data. Extending the efficient methods of 
statistical spectroscopy 11271, 12811 we exactly calculate the 
first and second moments of the Hamiltonian for differ- 
ent configurations at fixed spin and parity. As a practical 
tool we use the exact decomposition of many-body con- 
figurational space into classes corresponding to different 
parity and number of harmonic oscillator excitations. 
An accurate estimate of the shell-model ground state en- 
ergy is required being generally as time consuming as 
the previous steps. This stage can be improved by using 
the exponential convergence method suggested and ap- 
plied in Refs. ll^l^Bltl . or/and the recently developed 
projected configuration interaction method 

MM- In 

reverse, some knowledge about the level density can be 
helpful for extracting the ground state energy. 

The code described in this paper is based on nuclear 
statistical spectroscopy l28ll . It allows one to calcu- 
late the spin- and parity-projected moments of the nu- 
clear shell-model Hamiltonian, which can be further 
used for an accurate description of the level density up 
to about 15 MeV excitation energy. It can be also ap- 
plied to other mesoscopic systems, such as interacting 
cold atoms in harmonic oscillator traps. The code is 
para llelized using the Message Passing Interface (MPI) 
1 34J] and a master-slaves dynamical load-balancing ap- 
proach. The parallel code was thoroughly tested on 
the massively parallel computers at NERSC [35], and it 
shows very good scaling when using up to 4000 cores. 

The paper is organized as follows. In Sec. II the 
method of fixed spin- and parity-dependent configura- 
tion moments is revisited. The method allows to trace 
such quantum numbers as parity and angular momen- 
tum explicitly. The extension of the algorithm to the 
proton-neutron formalism is discussed in Sec. III. In 
Sec. IV we introduce the structure of the program and 
supply the examples of input files. Examples of calcu- 
lations are presented and compared to exact shell model 
results in Sec. V. Section VI is devoted to conclusions. 



2. Theory outline 

In this work we closely follow the approach proposed 
in Ref |25] (see also Refs. |22, 23t]). For clarity we 
repeat here the main ideas and equations we are going 
to use for calculating the level density. 

According to the method of moments one can cal- 
culate the density p(E, a) of levels with a given set of 
quantum numbers as a function of excitation energy E 



p(E, a) = ^ D aK ■ G aK (E), 



(D 



where a — {«, J, T z , n) includes all quantum numbers of 
interest, namely the number of particles n (protons and 
neutrons), total spin J, isospin projection T z , and par- 
ity 7r. The subscript k represents a configuration of n 
particles distributed over q spherical single-particle or- 
bitals. Each configuration k is fixed by a set of occupa- 
tion numbers, k — [k\,K2,... , K q }, where kj is the number 
of particles occupying the spherical single-particle level 
j. The configuration has a certain number of particles, 
total isospin projection, and parity. The sum in Eq. ([1} 
runs over all possible configurations corresponding to 
given values of n, T z , and n. The dimension D aK equals 
the number of correctly antisymmetrized many-fermion 
states with given J that can be built for a given config- 
uration k. The function G aK is a finite-range Gaussian 
defined as in 1E2I1 : 

G aK {E) - G(E + £g. s . - E aK , cr aK ), (2) 
exp(-;t 2 /2cr 2 ) , \x\<t}-<t 



G(x, a) = N ■ 







>TJ- 



(3) 



where E aK and cr aK are the fixed-/ centroids and widths, 
which will be defined later, E g . s . is the ground state en- 
ergy, 77 is the cut-off parameter, and N is the normal- 
ization factor corresponding to the following condition: 
J x G(x, cr)dx - 1 . 

A very important ingredient of the method is the ac- 
curate knowledge of the ground state energy £g. s .. It 
is also necessary to find an optimal value of the cut-off 
parameter rj, see the discussion in [25]. 

Assuming a two-body shell-model Hamiltonian, 



(4) 



ijkl 



we have to calculate traces of the first and second power 
of this Hamiltonian, Tr[H] and Tr[H 2 ], for each config- 
uration k which will determine the fixed- J centroids and 
widths in Eq. ©: 



where 



(H) a 
( Hl ) 



--Tv^[H]/D a , 
Tr iOK) [H 2 ]/D a 



(5) 
(6) 

(7) 
(8) 



If the many-body states |v, J) with a certain set of quan- 
tum numbers a = {v, J], including spin J, form a com- 
plete set for the configuration a:, the symbol of trace, 
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where i, j, I, and q are single-particle states with certain 
spin projections and possible occupation numbers equal 
Oorl. 

Notice that the single-particle orbitals we have used 
to define the configurations in Eq. ([TJ can host all par- 
ticles with all possible spin projections corresponding 
to spin of the orbital. The dimension factor D^ K = 
Tr^falaj] can be interpreted as a number of many- 
body states, possible for the configuration k, with the 
fixed projection M z (if we consider M z -traces) and un- 
der the condition that the single-particle state i is occu- 
pied: 

D™ = Tr«[aJaW], 



Figure 1: Speedup is defined as T\/T n , where T n is the calculation 
time, when n processors were used. These calculations were per- 
formed on the FRANKLIN supercomputer at the National Energy Re- 
search Scientific Computing Center (NERSC) (Hi. 



Tr m '[' • •], means the sum of all diagonal matrix ele- 
ments, 2 ( V ,J\ ' ' ' \v, /), within this subspace. Techni- 
cally, it is more convenient to derive these traces in a ba- 
sis with a fixed spin projection |v, M z ), Tr (M;) [- • •], rather 
than in the basis with fixed total spin |v, J), Tr (y) [- • •]. 
/-traces can be easily expressed through the M.-traces, 
given the rotational symmetry of the Hamiltonian, 

Tr (y) [---] = Tr w [••■]! -Tr ( ^ [•••]! . (9) 

'M Z =J <M z =J+\ 

For simplicity, in Eq. © we omitted all quantum num- 
bers, except the projection M z and the total spin /. 

Hereafter we use the label a to denote a set of quan- 
tum numbers that includes either the fixed M z or the 
fixed J, keeping in mind that Eq. (01 can always con- 
nect them. In every important case we will point out 
which set of quantum numbers is used. Following the 
approach of 13611 . we can obtain the following expres- 
sions for the traces in Eqs. ((7} and 

Tr^[#] = Yj * D « + Z V W D ™ ' (10) 

i i<i 

Tr M [ff J ] = Xl e ' ,D ™ + 



E 



q<l 



D m + 



{i<mq<D 



D m + 
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D [ tJ K q] = W aK >[a\aW q a q a jai \, 



Jok) 



(12) 



These D-structures were called propagation functions in 
i36Tl . For completeness, we repeat here the recipe used 
for calculating them. One can show [36] that 



D ln> 



• rJ = Z(-^ z 



(13) 



where all f, are non-negative integers, the configuration 
k' = {/c'j , k' v ... , K q } can be derived from the original con- 
figuration k - {k\,K2, ... ,K q ) by removing t particles 
from the single-particle states r\, r2, ■ ■ ■ r s . A formal ex- 
pression for the new configuration k! can be written as 
follows: 



k 'j = k j- z 



ti, 



(14) 



i (nej) 



where the sum includes only those values of i for which 
the corresponding single -particle state r, belongs to the 
single-particle level j. We also assume that all the oc- 
cupation numbers id- must be positive, which imposes 
certain restrictions on the possible values of the ampli- 
tudes ti. For every new configuration k' one can easily 
define new quantum numbers, a' = {«' 'M' z T'„tt'}, enter- 
ing Eq. ( fT3l . Examples are the new number of particles 
n' = n — t and the new spin projection, 



M' z = M z — t\tn,- x — t2tn r2 — ■ ■ ■ — t s m rs , 



(15) 



where m n is the M z projection of the single-particle state 
r,. The new isospin T' and parity n' are defined simi- 
larly. 

3. The algorithm for the method of moments in the 
proton-neutron formalism 

Here we describe some technical features of the algo- 
rithm developed for the calculation of the nuclear level 



density. We treat protons and neutrons separately, so 
that the basis of many-body wave functions is repre- 
sented by a product of proton and neutron parts: 

\v,M z ) = \v p ,M? ) )-\v n ,M™), (16) 

where Mi p) + M 7 n) = M 7 . The wave functions (IT6b have 
the fixed isospin projection T z , but do not have a certain 
isospin T. As we already mentioned, it is more con- 
venient to use the basis of the wave functions with the 
fixed spin projection M z , rather than the one with the 
fixed spin J. 

One could gain essential advantages from such a sep- 
aration of the basis. One of them is related to the num- 
ber of configurations that appear in the sum of Eq. ([1}. 
Naturally, the number of configurations with fixed T z 
is much greater than the number of configurations with 
fixed isospin. This allows the use of many-cores com- 
puters with greater efficiency. In other words, the cal- 
culation of the sum in Eq. ([TJ with a larger number of 
configurations can be more efficiently distributed over 
a larger number of processors. Fig. [T] presents the 
speedup (calculation speed gain) as a function of the 
number of used processors. One can see that the case 
with the larger number of configurations, 68 Se, scales 
better than the case with the lower number of configu- 
rations, ^Ge. Up to 2000 cores, the speedup is almost 
perfect (the dotted line presents an ideal speedup). At 
this point the calculation time is about 1-2 minutes and 
further improvement is hardly achievable. 

Another significant advantage of the proton-neutron 
formalism is the new algorithm of calculating the di- 
mensions D aK , D^l, Da{r, etc. Because of the proton- 
neutron separation one can calculate all proton and neu- 
tron dimensions separately. Later, the dimensions we 
are interested in can be easily constructed from the pro- 
ton and neutron parts using the convolution, 

M^ ) +M i " ) =M- 

where, instead of the whole set of quantum numbers 
a, only the spin projection M z was explicitly indicated. 
Here k p and k„ are the proton and neutron parts of the 
configuration k. Eq. dTvT > can be easily applied to all 
types of dimensions, D a '~}, needed in the formalism of 
section II. The advantage comes from the fact that one 
can calculate and keep in memory all proton and neutron 
dimensions, D M ( P ) K and D M w K , for all possible projec- 
tions mf and M ( "\ and for all possible configurations 
k p and k„. Afterwards, using Eqs. (fTTt and (TTTb . one can 
calculate very fast all the dimensions: D aK , D^l, D^J}, 
etc... , for all M z and J. 



One more technical detail, which allows a signifi- 
cant speed up of the algorithm, is that by using the 
proton-neutron separation one can avoid multiple com- 
putations of the most time consuming structures, such 
as Dal 9 ®. Let us consider a case when all four single- 
particle states {ijql} are protons. One can then use an 
equation similar to Eq. (fTTl i. 

= X D ^l/ D ^- (18) 

For all configurations k that have the same proton parts 
k b one would have to recalculate Z)''^' 1 for each neu- 

tron configuration. Alternatively, one can calculate 
0^5? on ly once, and store the results in memory. That 

M- Kp 

strategy, however, would require a large amount of stor- 
age. More efficiently, one can only store the contribu- 
tions of the D„i structures to the width, Eq. (fTTT i. that 
is, one can only store the following T-structures, 

(f<j)*(q<I) 

-4V qiil V qjfl ]D l ^l ; (19) 

where all single-particle states are protons. Thus, in- 
stead of using Eq. ([T8l > one can calculate the contribu- 
tion to the width directly via the convolution, 

Tr^[/7 2 ] = ...+ J] 7>V D m<»V ( 20 ) 

which is very similar to Eqs. ([PTt and ([T"8V The new ap- 
proach avoids multiple calculations of D jj?? ■ Storing 

the structures Eq. I ll 9) , one may significantly speed up 
the algorithm for large cases, such as 68 Se in pf + gg/2 
model space. The downside is that the calculation of the 
T-structures, 7 ' M w> K , T M (»> K , does not always scale well 
on a large number of cores, since the number of these 
T-structures is much smaller than the total number of 
configurations. 

4. Description of the program 

The program consists of two separate codes. The first 
code is called MM, which is the main code in the pro- 
gram. MM code performs calculation of the first and 
second moments for all the configurations within the 
given range of spins and for certain parity. It is the 
most complicated and resource demanding part of the 
program. MM requires parallel computing. The second 
code is very simple and fast. It takes the output of the 
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lar matrix elements, for example, to calculate Clebsch- 
Gordan coefficients; qsort.f contains quick-sort subrou- 
tines. 

The flowchart of the mm_1.15.f90 is shown in Fig|2l 
and some important subroutines in the mm_1.15.f90 are 
listed below. 

The subroutine prepare -interaction calculates the 
T-structures according to Eq. (TT9V This is an important 
part of the code which allows to speed up the program 
significantly. 

The subroutine cc_density_calc contains the main 
loop over the configurations (the loop back in Fig. 2), 
calculates the first and second moments based on Eqs. 
(17181101201 1. and saves the results. This subroutine uses 
the J-structures precalculated in the subroutine pre- 
pare -interaction. 

Both subroutines require parallel computing. 
The simplest "Master-Slave" parallel programming 
paradigm with Dynamic Load Balancing was used. 

4.2. MM input 

The input files include input.dat and the files that de- 
fine the single-particle model space as well as the inter- 
action in this model space. A typical input.dat file that 
specifies the parameters in the code is listed below. The 
parameters are followed by their meaning: 



Figure 2: The flowchart of the program. See text for details. 



MM code (the first and second moments) and builds the 
nuclear level densities according to Eqs. (Q]|3). It does 
not require parallel computing. 

Below we will concentrate only on description of the 
MM code. The detailed instructions on how to use the 
second code can be found in the readme.txt file, which 
is in the main project directory. 

4.1. The structure of the MM code 

In the MM code, the calculation of the first and sec- 
ond moments for a given nucleus is carried out. To 
compile this code simply follow the instructions in the 
readme.txt file (or type make) and then run the exe- 
cutable file mm.out. 

The code contains four files: mm_1.15.f90 (1.15 is 
the current version) contains the main subroutines and 
calls subroutines from other files to perform the calcu- 
lation; interaction.f90 reads the interaction file; angu- 
larme.for contains subroutines to work with the angu- 



66 IZN 

1 22 2 ! P 1J m i n 2y max 2y s t e p 

int/sd.spl ! single-particle model space 

int/usd.int ! interaction 

Here Z and N are the number of protons and neutrons, 
respectively, in the valence space; P is parity (P — 1 
corresponds to positive parity and P = -1 to negative 
parity); 2J m [ n , 2/ max , and 2J step define the range of to- 
tal spin for which the moments are to be calculated: the 
spin changes from the minimum value J m i n to the maxi- 
mum value 7 max with the step 7 step . In the example listed 
above the total spin changes from to 1 1 with the step 1 . 
The single-particle model space and the interaction are 
defined by two separate files and input.dat must have 
the names of these two files similar to the shown exam- 
ple. Detailed description of *.spl and *.int files is given 
in the readme.txt file. 

The above example describes the 28 Si nucleus in the 
sd-model space with the USD interaction. The mo- 
ments will be calculated for positive parity and for all 
possible spins from to 11. 
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4.3. MM output 

The main output of the MM code is presented in files 
mm_res_#.dat. These files are enumerated by spin num- 
ber #, for example mm_res_0.dat corresponds to spin 
7 = 0, mm_res_l.dat corresponds to spin J — 1/2, and 
so on (for more details see the readme.txt file). Each 
output file contains the data needed for the density cal- 
culation: dimensions, first and second moments. 

Another output file mm_conf.dat contains informa- 
tion about the configurations. This information is not 
used for the density calculation, but could be useful for 
checking and testing purposes. 

4.4. Density calculation 

After all the moments are prepared, the density can be 
calculated with the code den.out. The full description 
of density calculation can be found in the readme.txt 
file. We just repeat once more that there are several ex- 
ternal parameters that need to be prepared before the 
calculation. These parameters cannot be defined within 
the method of moments, namely: r\ - the cut-off param- 
eter, E g s . - the ground state energy, and the energy in- 
terval (E m i n , E max ) for calculating the level density. For 
more details see the readme.txt file. 

5. Examples 



Element 


Space 


Total dim 


Elapsed time (sec) 


70 Br 


Pf + 89/2 


10 15 


1.07- 10 4 


68 Se 


Pf + 89/2 


10 15 


1.03 • 10 4 


64 Ge 


Pf + 89/2 


10 14 


0.76 ■ 10 4 


60 Zn 


Pf 


10 11 


37.4 


52 Fe 


Pf 


10 10 


13.6 


28 Si 


sd 


10 6 


0.7 



Table 1 : Elapsed times of nuclear level density calculations (for all 
/, positive parity) with the moments method code. The calculations 
were done on a 16 cores machine with 2.8 GHz CPU frequency. 



Table[T]presents calculation times for different nuclei 
calculated in different shell-model spaces. The calcula- 
tions were done on a 16 cores machine with 2.8 GHz 
CPU frequency. One core ('master') distributed all the 
work between other 15 cores ('slaves'). One can em- 
phasize here that the listed times correspond to calcula- 
tions of the nuclear level densities for all J and for pos- 
itive parity. For the case of 68 Se the largest m-scheme 
dimension is about 10 15 . For each J the m-scheme di- 
mensions vary from 10 12 to 10 14 , which makes direct 



diagonalization impossible. Using the moments method 
and our algorithm we are able to calculate the shapes of 
nuclear densities for 68 Se in less then three hours on a 
16 cores machine. For a number of processors reaching 
one thousand, it will take only few minutes to complete 
the calculation. 



5.7. 28 Si, sd-model space 

As a first example we consider the level density for 
the 28 Si nucleus in the 5£/-shell model space, where we 
use the USD interaction [37|. Fig- [3]presents the com- 
parison of the exact shell-model level densities for dif- 
ferent spins (solid lines) with those obtained with the 
moments methods (dashed lines). 

Eqs. (Q]) and (f2| require the knowledge of the ground 
state energy E g . s . and the cut-off parameter r\. While 
the ground state energy of 28 Si can be calculated in this 
case using the standard shell model, £g S = -135.94 
MeV, for the value of the cut-off parameter 77 we have 
only a general idea that it should be around 3 II22[ 12311 . 
For a better description of level densities in the moments 
method we can adjust the 77 parameter to optimally re- 
produce the exact shell-model densities. From Fig. [3] 
one can see that choosing 77 = 2.8, the level densities 
of the moments method reproduce quite well the ex- 
act shell-model level densities. The cut-off parameter 
plays a role similar to that of the width in a Gaussian 
distribution. Indeed, if we increase the cut-off param- 
eter, the density becomes wider and lower, while de- 
creasing it leads to a narrowing of the density. One 
should also mention that the exact spin- and parity- 
dependent shell-model densities were calculated with 
the NuShellX code Q. 



5.2. M Ge, pf- and pf + g9 12— model spaces 

As mentioned in the Introduction, one could envi- 
sion using information from the level densities to extract 
with a good approximation the ground state energies. 
Using our algorithm and the moments method one can 
easily calculate the nuclear level density for any nucleus 
that can be described in the pf + gg/2 model space. The 
Hamiltonian used for this model space was built starting 
with the GXPF1 A interaction for the pf model space, to 
which the G-matrix elements that describe the interac- 
tion between the pf orbits and #9/2 orbit were added. 
The single-particle energy for the gg/2 orbit was fixed 
at -0.637 MeV. Even for the worst case, the calcula- 
tion takes about three hours for sixteen processors and 
only few minutes for one thousand processors. Fig. |4] 
presents the results obtained for 64 Ge, nucleus that is be- 
lieved to be a "waiting-point" along the rp-process path 
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Figure 3: Si, parity=+l. Comparison of nuclear level densities between the exact shell model (solid line) and the moments method (dashed line). 
Cut-off parameter rj = 2.8, interaction: USD, .sd-shell. 
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Figure 4: Level densities for the 64 Ge nucleus, J = 0, 2 and positive parity. Solid line presents the calculation in the p/-shell with the GXPF1A 
interaction. For this calculation we know the ground state energy, E gs _(pf) = -304.25 MeV. Other three lines present calculations in the large 
model space, when level gg/2 is added. The ground state energy for these cases is E gs _(pf + gg/2) = E gs {pf) - AE, where AE is the energy shift. 
The cut-off parameter is 77 = 2.6. 
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[39-41]. We only present the densities for J — 0, 2 and 
positive parity. 

The corresponding input.dat file for 64 Ge in 
pf -model space looks like 

12 12 \ZN 

1 22 2 ! P 2 J m [ n 2 J max 2 J s t e p 

int/pf.spl ! single-particle model space 

int/gxla.int ! interaction 

Here we have twelve protons and twelve neutrons in 
the pf -model space. The calculation is done for all 
spins between 7 = and 7=11 and positive parity. The 
single-particle space is defined in the pf.spl file and the 
interaction is in gxla.int. For the pf+g9/2 model space 
we need to change the single-particle file to pfg9.spl and 
the interaction file to pfg9.int. 

It is important to notice that in the pf model space 
the shell-model calculations of the ground state energies 
can be done. For 64 Ge in the pf '-shell we obtain the 
following ground state energy: 



E g . s .(pf) = -304.25 MeV. 



(21) 



Using this ground state energy and the cut-off parame- 
ter j] = 2.6, we are able to calculate the level densities 
according to Eqs. (Q} and (0. The solid lines in Fig. [4] 
represent the density in the pf -shell. 

To calculate the same level density in the pf + gg/2 
model space we have to adjust the ground state energy 
and the cut-off parameter for this space. For the cut-off 
parameter we use the same value, r\ = 2.6, but it is prac- 
tically impossible to calculate by shell-model diagonal- 
ization the ground state energy since the dimension is 
too large. The ground state energy for the larger model 
space, that is pf + gg/2, must be lower compared to the 
ground state energy for the smaller model space, that is 
p f. Let us introduce this energy difference, AE, as 



E g ,.(pf + ?9/ 2 ) = E g Apf) ~ AE. 



(22) 



The dotted lines in Fig. |4]show the level densities if we 
keep the ground state energy for pf + gg/2 model space 
as it was in the pf case keeping AE = 0. It is natural 
to expect only small differences between the level den- 
sities calculated in those two model spaces at low exci- 
tation energy since in the pf + gg/2 model space we use 
the same GXPF1 A interaction for the pf subspace. By 
decreasing the ground state energies for the pf + gg/2 
model space (introducing non-zero AE), one gets the 
dashed lines on Fig. [4] The dash-dotted lines there cor- 
respond to ground state energy E g . s .(j>/+ qg/2) = -305.8 
MeV of 64 Ge, which was obtained by a truncated shell 



model calculation with up to 6 particles excited from 
the fi/2 orbits and/or into the gg/2 orbit. The m-scheme 
dimension in this calculation, 13.5 x 10 9 , is at the up- 
per limit of the state of the art shell-model calculation. 
As one can see, this value does not describe satisfactory 
the level densities at low excitation energy. In order to 
make the low-lying part of the two densities very close 
(dashed and solid lines on Fig. 0J, one has to adjust the 
ground state energy for the pf + gg/2 model space to the 
following value: 



EgApf + 89,2) = -306.7 MeV for 



64 



Ge, (23) 



The "low-lying part of the density" should be chosen 
such that the excitations to the gg/2 orbit do not give 
a significant contribution. For these cases we use the 
interval 3-6 MeV in excitation energy. We conclude that 
the adjustment of Eq. ( 1231 can be treated as a method 
for estimating the ground state energies in larger spaces; 
for more details see D25J. 

6. Summary 

In summary, we have developed an efficient Fortran 
code for calculating the centroids and widths of the 
shell-model spin- and parity-dependent configurations, 
which can be used for calculating the nuclear level den- 
sities. The code is parallelized using the Message Pass- 
ing Interface (MPI) [34] and a master-slaves dynamical 
load-balancing approach. The parallel code was thor- 
oughly tested on the massively parallel computers at 
NERSC and it shows very good scaling when us- 
ing up to 4000 cores. The algorithm used takes advan- 
tage of the separation of the model space in neutron and 
proton subspaces. This separation provides two impor- 
tant advantages: (i) the exponentially exploding dimen- 
sions and propagators can be calculated more efficiently 
in proton and neutron subspaces, and the full results can 
be recovered via simple convolutions; (ii) the number of 
configurations is significantly increased in the proton- 
neutron formalism, considerably improving the scala- 
bility of the algorithm on massively parallel computers. 
Our tests indicate almost perfect scaling for up to 4000 
cores. The new algorithm is so fast that the bottleneck 
of the calculation is now that of the ground state en- 
ergy. That is why we could not test our algorithm for 
cases that take more than one minute on 4000 cores. 
Therefore, we investigated the possibility of using the 
calculated shapes of the nuclear level densities to ex- 
tract the ground state energy. We showed that by incre- 
menting the model space and the effective interaction, 
and imposing the condition that the level density does 
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not change at low expectation energy, one can reliably 
predict the ground state energy, and further the full level 
density. This new method of extracting the shell model 
ground state energy for model spaces whose dimensions 
are unmanageable for direct diagonalization opens new 
opportunities for calculating shell model level densities 
of heavier nuclei of interest for nuclear astrophysics, nu- 
clear energy and medical physics applications. 

A further development of the application of statisti- 
cal spectroscopy to nuclear level density is the removal 
of center-of-mass spurious states from the level density 
for shell model spaces that allow complete factorization 
of the center-of-mass and intrinsic wave functions. A 
new algorithm implementing this idea was recently pre- 
sented [42], and a high performance code was devel- 
oped. This new code could be made available upon re- 
quest. 
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